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Abstract 


Phase-plane  and  small-parameter  asymptotic  techniques  are  used  to  analyze  systems 
of  ordinary  differential  equations  that  describe  the  transient  behavior  of  non-Newtonian 
fluids  in  shear  flow.  These  systems  approximate  the  partial  differential  equations  that 
derive  from  three-dimensional  balance  laws  and  from  differential  constitutive  models  for 
highly  elastic  liquids.  Two  models  are  considered:  one  with  a  single  relaxation  time  and 
small  Newtonian  viscosity:  the  other  with  two  relaxation  times  and  no  Newtonian  viscosity. 
Both  possess  the  key  feature  that  the  variation  of  steady  shear  stress  with  strain  rate  is 
not  monotone. 


The  analysis  shows  that  both  models  exhibit  several  distinctive  phenomena:  spurt, 
shape  memory,  hysteresis,  latency,  and  normal  stress  oscillations.  The  predictions  for  the 
spurt  phenomenon  agree  quantitatively  with  experimental  results  for  polymer  melts:  the 
other  new  phenomena,  which  were  discovered  recently  in  numerical  simulation,  should  also 


be  observable  in  rheological  experiments. 
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1.  Introduction 

The  purpose  of  this  paper  is  to  analyze  novel  phenomena  in  dynamic  shearing  flows  of 
non-Newtonian  fluids  that  are  important  in  polymer  processing  [15].  One  striking  effect, 
called  the  “spurt”  phenomenon,  was  observed  by  Vinogradov  et  aj.  [17]  in  experiments  con¬ 
cerning  quasi-static  flow  of  monodispersive  polvisoprenes  through  capillaries.  They  found 
that  the  volumetric  flow  rate  increases  dramatically  at  a  critical  stress  that  is  independent 
of  molecular  weight.  Until  recently,  spurt  has  been  associated  with  failure  of  the  polymer 
to  adhere  to  the  wall  [5].  The  focus  of  our  current  research  is  an  alternative  explanation 
of  spurt  and  related  phenomena. 

The  present  paper  analyzes  systems  of  ordinary  differential  equations  that  approxi¬ 
mate  the  dynamics  of  one-dimensional  shear  flow  of  highly  elastic  non-Newtonian  fluids 
at  low  Reynolds  number.  The  analysis  reveals  several  distinctive  phenomena  related  to 
spurt:  shape  memory,  hysteresis,  latency,  and  normal  stress  oscillations.  Furthermore, 
the  comparison  two  different  viscoelastic  fluid  models  shows  that  the  effect  of  a  second 
relaxation  time  is  similar  to  that  of  Newtonian  viscosity. 

The  essential  aspects  of  spurt  experiments  are  reflected  in  one-dimensional  shear  flow 
through  a  slit  die,  which  is  governed  by  fluid  balance  laws  and  constitutive  relations.  The 
key  feature  of  constitutive  relations  that  exhibit  spurt  is  a  non-monotonic  relation  between 
the  steady  shear  stress  and  strain  rate.  This  causes  jumps  in  the  steady  strain  rate  to  form 
when  the  driving  pressure  gradient  exceeds  a  critical  value;  such  jumps  correspond  to  the 
sudden  increase  in  volumetric  flow  rate  observed  in  the  experiments  of  Vinogradov  et  a 1. 

The  constitutive  models  we  study  derive  from  three-dimensional  differential  constitutive 
relations  with  multiple  relaxation  times  (based  on  work  of  Johnson  and  Segalman  [8]  and  of 
Oldroyd  [14]).  These  models  lead  to  systems  of  partied  differential  equations  that  are  evo¬ 
lutionary,  axe  globally  well-posed  (in  a  sense  described  below),  and  possess  discontinuous 
steady  states  that  explain  spurt  effects. 

The  Johnson-Segalman-Oldroyd  shear  flow  equations  were  investigated  using  numer¬ 
ical  simulation  in  Refs.  [9,  11].  Transient  flows  at  high  Deborah  (Weissenberg)  number 
and  low  Reynolds  number  exhibited  spurt,  shape  memory,  and  hysteresis.  Moreover,  other 
effects  were  evident,  such  as  latency,  normal  stress  oscillations,  and  molecular  weight  de¬ 
pendence  of  hysteresis.  A  principal  conclusion  is  that  satisfactory  explanation  of  the  spurt 
phenomenon  requires  understanding  the  full  dynamics  of  the  flow,  not  only  steady  flows. 
Furthermore  it  was  shown  [11]  that  the  equations  reduce  to  planar  systems  of  ordinary  dif¬ 
ferential  equations  in  the  flow  regime  appropriate  to  the  spurt  experiments.  The  objective 
of  the  present  paper  is  to  study  the  reduced  dynamical  systems  using  analytical  methods 
and  to  verify  the  dramatic  implications  of  numerical  simulation. 

Earlier  work  has  drawn  connections  between  spurt  phenomena  and  non-monotone 
constitutive  relations.  Hunter  and  Slemrod  [7]  studied  the  qualitative  behavior  of  discon¬ 
tinuous  steady  states  in  a  simple  one-dimensional  viscoelastic  model  that  exhibits  spurt, 
shape  memory,  and  hysteresis.  Their  model  contrasts  with  the  ones  studied  here,  however,  '* 
in  that  the  effects  observed  in  their  model  are  related  to  the  linear  instability  and  loss  of 
evolutionarity  of  the  flow  at  certain  Deborah  numbers.  t/ 

The  two  systems  of  partial  differential  equations  that  we  use  to  model  shear  flow  ^  ,. 
through  slit  dies  are  formulated  in  Sec.  2.  For  the  first  system  the  total  shear  stress  .riJ/0r 
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is  decomposed  into  a  polymer  contribution,  evolving  in  accordance  with  a  differential 
constitutive  relation  with  a  single  relaxation  time,  and  a  Newtonian  viscosity  contribution 
[system  (75(9)  in  Sec.  2j.  The  second  system  is  based  on  a  differential  constitutive  law 
with  two  widely  spaced  relaxation  times  but  no  Newtonian  viscosity  [system  (7502)  in. 
Sec.  2].  We  shall  show  that  the  dynamics  of  these  two  models  is  similar;  this  suggests  that 
the  qualitative  behavior  we  find  is  not  special  to  the  particular  viscoelastic  models  we  have 
analyzed. 

For  the  highly  elastic  and  very  viscous  polyisoprenes  used  in  the  spurt  experiment, 
the  ratio  a  of  Reynolds  number  to  Deborah  number  is  several  orders  of  magnitude  smaller 
than  the  ratio  e  of  Newtonian  viscosity  to  shear  viscosity.  This  implies  that  the  flow 
equations  (750)  are  approximated,  at  each  fixed  position  x  in  the  channel,  by  the  following 
dynamical  system  [11]: 


i  =  (Z  +  1)  (—■'J 


(1.1) 


Here  the  dot  denotes  the  time  derivative  and  T  =  fz,  with  /  being  the  driving  pressure 
gradient. 

Although  it  has  not  been  proved  mathematically  that  solutions  of  system  (1.1)  ap¬ 
proximate  solutions  of  (750),  this  conclusion  is  reinforced  by  some  analytical  results  as 
well  as  by  numerical  experiments.  Nohel,  Pego  and  Tzavaras  [13]  established  that  dis- 
continous  steady  states  are  stable  in  a  simple  model  in  which  the  polymer  shear  stress 
satisfies  a  single  differential  constitutive  relation;  this  model  can  be  arranged  to  have  the 
same  behavior  in  steady  shear  as  system  (750).  Research  in  progress  suggests  that  this 
is  true  also  of  the  full  systems  (750)  and  ( JSO2 ),  at  least  when  a  is  sufficiently  small. 

Analysis  of  system  (1.1)  is  carried  out  in  Sec.  3.  We  show  that  steady  states  of  sys¬ 
tem  (750),  some  of  which  axe  discontinuous,  correspond  to  families  of  critical  points  of 
the  dynamical  system.  We  deduce  the  local  characters  of  the  critical  ooints,  and  prove  that 
the  system  has  no  periodic  orbits  or  closed  separatrix  cycles.  More*,  er  we  determine  the 
global  dynamics  completely  through  a  phase-plane  analysis;  this  is  a  •  nplished  by  con¬ 
structing  a  natural  Lvapunov-like  function.  Thus  we  identify  the  globally  asymptotically 
stable  steady  states  for  each  position  x. 

The  same  approximation  in  the  model  ( JSO2 )  leads  to  the  study  a  similar  but  more 
complicated  dynamical  system.  In  Sec.  4  we  show  that,  in  the  relevant  range  of  parameters, 
the  essential  dynamical  features  of  the  latter  system  are  qualitatively  the  same  as  those  of 
system  (1.1). 

In  Sec.  5  the  results  of  Secs.  3  and  4  is  applied  to  explain  spurt,  shape  memory,  hystere¬ 
sis,  and  other  effects  observed  in  the  numerical  simulations.  This  involves  understanding 
how  the  phase  portraits  of  system  (1.1)  change  as  T  is  varied.  In  Sec.  6  we  develop  asymp¬ 
totic  expansions  in  powers  of  e  of  solutions  of  the  dynamical  systems;  this  enables  us  to 
explain  the  latency  effect  (a  pseudo-steady  state  that  precedes  spurt).  The  asymptotic 
analysis  also  permits  a  more  quantitative  comparison  of  the  dynamics  of  the  two  constitu¬ 
tive  models  when  e  is  sufficiently  small.  In  Sec.  7  we  discuss  physical  implications  of  our 
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results,  particularly  those  predictions  that  warrant  testing  by  rheological  experiment. 

Finally,  in  Sec.  8,  we  draw  certain  conclusions:  although  the  analysis  in  this  paper 
applies  only  to  the  special  constitutive  models  we  have  studied,  we  expect  that  the  qual¬ 
itative  features  of  our  results  appear  in  a  broad  class  of  non-Newtonian  fluids.  Indeed, 
numerical  simulation  by  Kolkka  and  Ierley  [10]  shows  that  another  model  with  a  single 
relaxation  time  and  Newtonian  viscosity  exhibits  similar  character. 

2.  Mathematical  Formulation  of  Generalized  Johnson-Segalman  Models 

The  motion  of  a  fluid  under  incompressible  and  isothermal  conditions  is  governed  by 
the  balance  of  linear  momentum 


=  VS.  (2.1) 

Here  p  is  the  fluid  density,  v  is  the  particle  velocity,  and  S  is  the  total  stress  tensor.  The 
response  characteristics  of  the  fluid  are  embodied  in  the  constitutive  relation  for  the  stress. 
For  viscoelastic  fluids  with  fading  memory,  these  relations  specify  the  stress  as  a  functional 
of  the  deformation  history  of  the  fluid.  Many  sophisticated  constitutive  models  have  been 
devised:  see  Ref.  [2]  for  a  survey.  In  the  present  work,  we  focus  on  a  class  of  multi-mode 
differential  models.  These  models  generalize  the  Johnson-Segalman  model  [8]  by  allowing 
for  a  distinct  slip  parameter  for  each  mode;  alternatively,  one  can  regard  each  mode  to 
be  governed  by  a  special  case  of  the  Oldroyd  constitutive  equation  [14].  As  discussed  in 
Secs.  1  and  8,  we  believe  that  qualitative  aspects  of  our  results  are  not  limited  to  these 
particular  models. 

To  specify  these  constitutive  relations,  we  decompose  the  total  stress  as 

S  =  -pi  -f  2t?D  +  n  .  (2.2) 

In  this  equation,  p  is  an  isotropic  pressure  (which  is  determined  from  the  incompressibility 
constraint)  and  q  is  the  coefficient  of  Newtonian  viscosity;  we  let  D  :=  I  [Vv  -f  (Vv  )r] 
and  Q  :=  \  [Vv  —  (Vv  )r]  be  the  symmetric  and  antisymmetric  parts  of  the  velocity  gra¬ 
dient  Vv,  which  has  components  (Vv  ),J  :=  dv'  jdx^ .  The  tensor  II  is  the  non-Newtonian 
extra  stress,  which  we  decompose  into  partial  stresses: 


dv 

dt 


+  v  •  Vv 


n  :=  £  • 

1=  1 


(2.3) 


Each  partial  stress  t,  is  specified  by  the  differential  constitutive  law 


r;  —  2|U;D  —  A  jXj  , 


(2.4) 


where 


dx. 


+  (v  •  V)xj  -  XjVv  -  (VvfjTj  +  (1  -  aj)(xjD  +  D xj)  (2.5) 


-4- 


is  the  objective  time  derivative  of  xj  with  slip  parameter  a;.  The  parameters  are  elastic 
shear  moduli,  and  the  A;  axe  relaxation  rates;  the  quotients  Pj/\j  axe  shear  viscosities. 

Essential  properties  of  constitutive  relations  are  exhibited  in  simple  planar  Poiseuille 
shear  flow.  With  the  flow  aligned  along  the  y-axis  (see  Fig.  1),  the  flow  variables  are 
independent  of  y.  Therefore  the  velocity  field  is  v  =  (0,u(x,t)),  and  the  balance  of  mass 
is  automatically  satisfied.  Furthermore,  the  components  of  each  partial  stress  tensor  xj 
can  be  written  tj1  =  7 xxy  =  x'j1  =  <7j(x,f),  and  =  Tj(x,f),  while  the  pressure 
takes  the  form  p  =  po(x,  t)  —  f(t)y,  f  being  the  pressure  gradient  driving  the  flow.  In  these 
terms.  Eqs.  (2.4)  become 

Ijt  +  (1  -  aj)<Tjvx  =  - \jlj  , 

ajt  ~  [y(l  +  o.j )lj  -  y(l  -  Q;)tj  +  V-i]  v *  =  >  (2-6) 

Tjt  -  (I  +  a.j)cTjVx  =  -A jTj  . 


Fig.  1:  Shear  flow  through  a  slit-die. 

For  each  j ,  we  introduce  the  variables 

zj  :=  y(!  +  ai )lj  ~  i(x  ~  Qj)rj  (2-7) 
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and 

Wi  :=  ~l(l  +  ai hi  ~  l(X  ~  ai)Ti  •  (2-3) 

Taking  linear  combinations  of  the  first  and  last  of  Eqs.  (2.6)  yields  an  equivalent  pair 
of  equations  for  Zj  and  Wj.  The  equation  for  Wj  is  Wjt  =  —\jWj\  noting  that  Wj 
must  remain  finite  a s  t  — oo,  we  conclude  that  Wj  =  0.  Thus  the  two  evolution 
equations  for  7 j  and  Tj  can  be  replaced  by  a  single  equation  for  Zj.  Furthermore,  Zj  — 
—  i(l  —  aj)(Tj  —  7 j),  so  that  Zj  represents  a  contribution  to  the  principal  normal  stress 

difference,  j^jLl  (xjy  -  xfx)  =  ( Tj  ~  7 j)- 

With  these  observations,  Eqs.  (2.6),  combined  with  the  balance  of  linear  momentum, 
Eq.  (2.1),  reduce  to 

m 

pvt  -  =  Vvxx  +  f  > 

j=1  (2.9) 

<7jt  -  (Zj  +  fij)  Vr  =  -A jffj  , 

Zjt  +  (l  -  a])  crjvx  =  -\jZj  . 

This  system  serves  as  a  basis  for  the  flow  problems  discussed  below.  For  later  reference, 
notice  that  the  momentum  equation  in  Eq.  (2.9)  can  be  written  as  pvt  —  Tr  =  /,  where 

m 

T  :=  <Tj  +  t]vx  (2-10) 

i=i 

is  the  total  shear  stress.  If  77  =  0,  it  is  readily  shown  that  system  (2.9)  is  hyperbolic  if 
J2J ~i(Z}  +  pj)  >  0,  with  wave  speeds  ±  p~l  ^  0  (repeated  2m  -  1 

times).  Moreover,  for  initial-boundary-value  problems  appropriate  for  shear  flow,  and  for 
smooth  data  compatible  with  the  boundary  conditions,  techniques  in  Ref.  [16]  can  be 
used  to  establish  global  well  posedness  of  classical  solutions  if  the  data  have  small  total 
variation,  and  finite-time  blow-up  of  classical  solutions  if  the  data  have  large  variation. 

If  p  >  0,  general  theory  developed  in  Ref.  [13,  Sec.  3  and  App.  A]  yields  global 
existence  of  classical  solutions  for  smooth  initial  data  of  arbitrary  size,  and  also  existence 
of  almost  classical,  strong  solutions  with  discontinuities  in  the  initial  velocity  gradient  and 
in  the  stress  components;  the  latter  result  allows  one  to  prescribe  discontinuous  initial  data 
of  the  same  type  as  the  discontinuous  steady  states  studied  in  this  paper. 

In  what  follows,  we  will  consider  two  special  cases:  (A)  the  case  of  a  single  relaxation 
time  (m  =  1)  in  the  presence  of  Newtonian  viscosity  (77  ^  0);  and  (B)  the  case  of  two 
relaxation  times  (m  =  2)  but  no  Newtonian  viscosity  (77  =  0).  The  study  of  these  two 
models  determines  how  the  effect  of  Newtonian  viscosity  compares  to  the  effect  of  a  second 
relaxation  time.  In  both  cases  we  study  the  Poiseuille  shear  flow  between  parallel  plates 
located  at  x  —  ihf 2  (cf.  Fig.  1). 
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A.  Single  Relaxation  Time.  To  simplify  notation,  we  omit  the  subscript  j,  and  we 
nondimensionalize  the  variables  by  scaling  distance  by  h ,  time  by  A-1,  and  stress  by 
Furthermore,  if  we  replace  a,  v ,  and  /  by  a  (1  —  a2)l^2a,  v  :=  (1  —  a2)1//2u,  and 
/  :=  (1  —  a2)1/2/,  respectively,  then  the  parameter  a  disappears  from  Eqs.  (2.9).  Since  no 
confusion  will  arise,  we  omit  the  caret.  There  are  two  essential  dimensionless  parameters: 

a  :=  ph}  A2/**  ,  (2.11) 

a  ratio  of  Reynolds  number  to  Deborah  number;  and 

£  :=  rjX/fi  ,  (2.12) 


a  ratio  of  viscosities. 

The  resulting  initial-boundary- value  problem  governing  the  flow  [11]  is  the  system 

OcV i  (7 x  —  “1“  f  j 

crt  -  (Z  +  l)ur  =  ~<j  ,  ( JSO ) 

Zt  +  o"vx  =  —Z 

on  the  interval  [-1/2,0],  with  boundary  conditions 

u(-l/2,f)=0  and  yr(0,f)  =  0  (£C) 

and  initial  conditions 

u(x,0)  =  ”o(x)  ,  <r(x,0)  =  cr0(x)  ,  and  Z(x,  0)  =  Z0(x)  ,  (IC) 

where  y0(-l/2)  =  0,  Uq(0)  =  0  and  <7O(0)  =  0.  When  £  =  0,  and  Z  +  1  >  0,  the  system  is 
hyperbolic,  with  characteristics  speeds  0  and  ±[(Z  +  l)/a]1,/2  and  0. 

The  steady-state  solutions  of  system  {JSO),  when  the  forcing  term  /  is  a  constant 
/,  plav  an  important  role  in  our  discussion.  Such  a  solution,  denoted  by  u,  a,  and  Z,  can 
be  described  as  follows.  The  stress  components  a  and  Z  are  related  to  the  strain  rate  vx 


through 

_ 

a  =  - 7 

1  +  V' 

(2.13) 

and 

Z  +  1  —  .  • 

1  +  vx 

(2.14) 

Therefore,  the  steady  total  shear  stress  T  :=  cr  +  evx  is  given  by  T  — 

w(Ux),  where 

U'(j):-l  +  s2+“- 

(2.15) 

The  properties  of  w.  the  steady-state  relation  between  shear  stress  and  shear  strain 
rate,  are  crucial  to  the  behavior  of  the  flow.  By  symmetry,  it  suffices  to  consider  s  >  0. 
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For  all  £  >  0,  the  function  w  has  inflection  points  at  3  =  0  and  s  =  y/3.  When  e  >  1/8, 
the  function  w  is  strictly  increasing,  but  when  £  <  1/8,  the  function  w  is  not  monotone. 
Lack  of  monotonicity  is  the  fundamental  cause  of  the  non-Newtonian  behavior  studied  in 
this  paper,  so  hereafter  we  assume  that  £  <  1/8. 

The  graph  of  w  is  shown  in  Fig.  2.  Specifically,  w  has  a  maximum  at  s  =  s m  and  a 
minimum  at  s  =  sm,  where 


■®  M  i  m 


1  -  2s  T  Vl-Se 


ll/2 


(2.16) 


respectively,  at  which  points  it  takes  the  values  T m  :=  w(s\f)  and  Tm  •*=  w(sm).  The 
function  w  has  an  inflection  point  at  s-lnfl  =  \/3.  As  £  — ►  1/8  from  below,  the  two  extrema  of 
w  coalesce  at  Another  limit,  £  — ►  0,  occurs  in  applications  to  rheological  experiments 
(see  Sec.  6);  in  this  vein,  we  record  the  following  asymptotic  formulae: 


sm  =  1  2s  +  6s2  +  0(s3)  ,  (2-17) 

Tm  =  1/2  +  s  +  s2  +  0(s3)  ,  (2.18) 

=  e-'/'fl  -  (3/2)e  -  (2S/S)e2  +  0(e3)]  ,  (2.19) 

Tm  =  £1/2(2  -  £  -  (5/4)g2  +  Ofs'  jl  .  (2.20) 


The  momentum  equation,  together  with  the  boundary  condition  at  the  centerline, 
implies  chat  the  steady  total  shear  stress  satisfies  T  =  —fx  for  every  x  6  [—  0],  Therefore, 

the  steady  velocity  gradient  can  be  determined  as  a  function  of  x  by  solving 

v)(yx)  =  -lx  .  (2.21) 

Equivalently,  a  steady  state  solution  vx  satisfies  the  cubic  equation  P{vx)  —  0,  where 

P{s)  :=ss3-Ts2  +  (l  +  £)s-T.  (2.22) 


The  steady  velocity  profile  in  Fig.  3  is  obtained  by  integrating  vx  and  using  the  boundary 
condition  at  the  wall.  However,  because  the  function  iv  is  not  monotone,  there  might 
be  up  to  three  distinct  values  of  vx  that  satisfy  Eq.  (2.21)  for  any  particular  x  on  the 
interval  [—1/2,0].  Consequently,  vx  can  suffer  jump  discontinuities,  resulting  in  kinks  in 
the  velocity  profile  (as  at  the  point  xn  in  Fig.  3).  Indeed,  a  steady  solution  must  contain 
such  a  jump  if  the  total  stress  T wail  =  // 2  at  the  wall  exceeds  the  total  stress  Tm  at  the 
local  maximum  M  in  Fig.  2. 

Let  us  describe  in  more  detail  the  analytical  solutions  for  steady  flow  fields  in  two 
instances:  subcriticol  flow  (having  no  kinks)  and  supercritical  flow  with  one  kinkj_as  is 
pictured  in  Fig.  3).  To  this  end,  it  is  useful  to  define  s,np(T)  for  Tm  <T  and  T  <  — to 
be  the  root  s  of  w(s)  =  T  with  largest  absolute  value,  while  3Sub(T)  for  —T m  <T<  Tm 
is  the  root  with  smallest  absolute  value.  Then  we  define  the  following  antiderivatives: 


_ i  r!z 

^supC1)  ;==  /  J  s9np(T)  dT  , 

_ i  r~~fz 

V„ „(*):-/  J  5,«b(T)  dT  . 


(2.23) 
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w 


i 


Fig.  2:  Total  steady  shear  stress  T  vs.  shear  strain  rate  vx  for 
steady  flow.  The  case  of  three  critical  points  is  illustrated;  other 
possibilities  are  discussed  in  Sec.  3. 


Because  T(x )  =  — /x,  the  steady  velocity  field  is  given  as  follows: 


subcritical  flow: 
supercritical  flow: 


v(x)  =  Vsuh(x)-Vsah(-l/2)  ; 


v(x) 


K>up(x)  ~  v;up(— 1/2)  ,  X  <  x„  , 

lsup(^*)  “b  hsub(^)  i  X  X*  , 


(2.24) 

(2.25) 


where  x.  is  the  location  of  the  jump  in  vx,  which  is  related  to  the  stress  X,  at  the  discon¬ 
tinuity  by 

x.  =  -§T./Tw»u  .  (2.26) 

Determination  of  x.  and  X.  from  the  system  dynamics  is  discussed  in  Sec.  6.  Integrating 
the  velocity  with  respect  to  x  yields  the  flow  rate  per  cross-section 

Q  :=  2  f  v(x)  dx  .  (2.27) 

J- 1/2 
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Fig.  3:  Velocity  profile  for  steady  flow. 


This  quantity  is  the  one  most  easily  observed  experimentally  and,  along  with  Twau,  is  the 
primary  datum  describing  the  spurt  phenomenon  reported  by  Vinogradov  et  aJ.  [17]. 

This  approach  to  obtaining  analytical  steady  flows,  if  applied  directly  using  the  cubic 
formula,  yields  an  expression  for  v  that  is  intractable.  An  alternate  approach  is  to  integrate 
the  asymptotic  expression  for  vs  that  is  obtained  in  Sec.  6.  This  is  straight-forward  and 
results  in  a  fairly  compact  closed-form  formula  [18].  On  the  other  hand,  changing  variables 
of  integration  from  from  shear  stress  to  strain  rate  in  Eqs.  (2.23),  yields  integrals  that  can 
be  evaluated  explicitly;  we  refer  to  Ref.  [IS]  for  details  of  this  development. 

B.  Two  Relaxation  Times.  By  analogy  with  the  previous  model,  we  nondimensionalize 
the  variables  by  scaling  distance  by  h,  time  by  A]"1,  and  stress  by  ;  also,  we  replace  <t;, 
y,  and  /  by  a }  :=  (1  -  aj)1/'2<rJ1  v  (1  —  al)1/,2u,  and  /  :=  (1  -  af)1^2/,  respectively,  and 
we  omit  the  caret.  It  proves  convenient  to  choose  dimensionless  variables  as  follows: 


II 

x> 

3- 

V* 

»-*  KJ 

T: 

(2.28) 

=  ^2/^1 1 

(2.29) 

=  (l-a2)/(l-a2)  , 

(2.30) 

=  Ai^2/(A2Ml)  > 

(2.31) 
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In  this  notation,  the  initial-boundary-value  problem  for  the  system  becomes 

&Vt  —  [<Ji  +  c2]x  =  f  > 

Cl t  ~  {Z\  +  l)ur  =  ~al 

Zu  +  civx  =  -Zi,  (JSO2) 

<7 2t  ~  {Z2  "b  fl)vx  —  ~K  1 0C2  , 

Z2t  4"  8c2Vx  ~  tc  fiZ 2 

on  the  interval  [-1/2,0],  with  boundary  conditions  (BC)  and  initial  conditions  v(x,0)  = 
vq(x)  and 

<7j(x,  0)  =  coj(x)  and  Zj(x,0)  =  Zoj(x)  ,  i^j) 

where  j  =  1,2:  again  we  require  that  vq{—  1/2)  =  0,  Uq(0)  =  0,  and  cr0j(0)  =  0.  The 
condition  for  hvperbolicity  of  system  ( JS02 )  is  +  1  +  Z2  +  /?  >  0,  in  which  case  the 
characteristics  speeds  are  0  (repeated  three  times)  and  ±[(Z  1  +  1  +  Z2  +  0)/a]1/2. 

The  steady-state  solutions  of  system  ( JSO2 )  are  given  by 


c  1 

a2 


(2.32) 


Z 1+1  — 


Z2  +  0  — 


1  +  v: 


1  -|-  8rz2fi~2vl 


(2.33) 


Again,  the  steady  total  shear  stress  T  =  ax  -f  cr2  satisfies  T  —  —fx\  thus  the  steady  state 
solution  uz  satisfies  the  quartic  equation  P2(vx)  =  0,  where 

P2(s )  :=  -T5k2$~2s4  +  («  +  SK2i3~2)s:i  -  T(  1  -f  8k2(3~2)s 2  +  (1  +  k)s-T  .  (2.34) 


We  note  that  if  5  =  0  and  *  =  e,  then  P2  =  P,  and  the  steady  states  of  system  ( JS02 ) 
and  of  system  (JSO)  coincide. 
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3.  Phase  Plane  Analysis  for  System  ( JSO )  When  a  =  0 

A  great  deal  of  information  about  the  structure  of  solutions  of  system  (JSO)  can  be 
garnered  by  studying  a  system  of  ordinary  differential  equations  that  approximates  it  in 
a  certain  parameter  range.  Motivation  for  this  approximation  comes  from  the  following 
observation:  in  experiments  of  Vinogradov  et  a i.  [17],  a  =  ph2\2//j  is  of  the  order  10~12; 
thus  the  term  avt  in  the  momentum  equation  of  system  (JSO)  is  negligible  even  when  vt  is 
moderately  large.  We  are  led  to  study  the  approximation  to  system  (JSO)  obtained  when 
a  =  0.  The  behavior  of  solutions  of  the  resulting  dynamical  system  offers  an  explanation  for 
several  features  of  the  solutions  of  the  full  system  (JSO)  observed  in  the  computations  of 
Refs.  [9,  11];  in  fact,  these  calculations  prompted  the  following  analysis,  which  determines 
the  dynamics  of  the  approximating  system  completely. 

When  a  =  0,  the  momentum  equation  in  system  (JSO)  can.  be  integrated,  just  as  in 
the  case  of  steady  flows,  to  show  that  the  total  shear  stress  T  :=  a  +  evx  coincides  with 
the  steady  value  T(x)  =  —  fx.  Thus  T  =  T(x)  is  a  function  of  x  only,  even  though  cr  and 
vz  are  functions  of  both  x  and  t.  The  remaining  equations  of  system  (JSO)  yield,  for  each 
fixed  x,  the  autonomous,  quadratic,  planar  system  of  ordinary  differential  equations 


(3.1) 


Here  the  dot  denotes  the  derivative  d/dt.  We  emphasize  that  for  each  /,  a  different 
dynamical  system  is  obtained  at  each  x  on  the  interval  [—1/2,0]  in  the  channel  because 
T  =  —  fx.  By  symmetry,  we  may  focus  attention  on  the  case  T  >  0;  this  is  assumed 
throughout  Secs.  3  and  4.  The  dynamical  system  (3.1)  can  be  analyzed  completely  by  a 
phase-plane  analysis,  which  we  carry  out  in  detail.  Recall  from  Sec.  2  that  we  assume 
e  <  1/S  throughout. 

The  critical  points  of  system  (3.1)  satisfy  the  algebraic  system 


(Z  +  l  +  e)(=-l)+£-0, 


c  T\T  J 


These  equations  define,  respective!*'  a  hype  'bola  and  a  parabola  in  the  a-Z  plane;  these 
curves  are  drawn  in  Fig.  4.  which  >  -responds  to  the  most  comprehensive  case  of  three 
critical  points.  The  critical  point  intersections  of  these  curves.  In  particular,  critical 
points  lie  in  the  strip  0  <  cr  <  T. 

Eliminating  Z  in  these  equation.-  nows  that  the  cr-coordinates  of  the  critical  points 
satisfy  the  cubic  equation  Q(<r  / 1)  =  0,  where 


«(«):=  —«?-!)  +  !+«  «-!)+«• 
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Fig.  4:  The  phase  plane  in  the  case  of  three  critical  points. 

A  straightforward  calculation  using  Eq.  (2.22)  shows  that 

P(vr)  =  P  =  .  (3.4) 

Thus  each  critical  point  of  the  system  (3.1)  defines  a  steady-state  solution  of  system  ( JSO ): 

such  a  solution  corresponds  to  a  point  on  the  steady  total-stress  curve  (see  Fig.  2)  at  which 
the  total  stress  is  T(x).  Depending  on  the  value  of  T,  and  hence  x,  there  are  either  one,  two, 
or  three  critical  points.  To  simplify  the  ensuing  discussion,  we  shall  ignore  the  degenerate 
cases  of  two  critical  points,  where  T  —  Tm  or  T  =  Tm.  Consequently,  we  have: 

Proposition  3.1:  For  (almost )  all  positions  x  in  the  channel,  there  are  three  possibilities: 

(1 )  there  is  a  single  critical  point  A  when  T(x)  <  Tm; 

(2)  there  is  a  single  critical  point  C  when  T(x)  >  ; 

(3)  there  are  three  critical  points  A,  B,  and  C,  with  cr a  >  & B  >  &c,  when  Tm  <  T(x )  < 
Tm. 

To  determine  the  qualitative  structure  of  the  dynamical  system  (3.1),  we  first  study 
the  nature  of  the  critical  points.  The  behavior  of  orbits  near  a  critical  point  depends  on 
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the  linearization  of  Eq.  (3.1)  at  this  point,  i.e.,  on  the  eigenvalues  of  the  Jacobian 

(3.5) 

evaluated  at  the  critical  point.  The  character  of  the  eigenvalues  of  J  can  be  determined 
from  the  signs  of  the  trace  of  J,  given  by 

— £  Tr  J  =  Z  +  1  +  2s  ;  (3 .6) 

the  determinant  of  J,  given  by 

eDetJ  =  Z  +  l  +  £  +  E^2=-l)  (=-l)  ; 

and  the  discriminant  of  J,  given  by 

£2  Discrm  J  =  (Z  +  l)2  -  8 T2  (=  -  +  §T2  . 

We  note  a  useful  fact:  at  a  critical  point, 

eDetJ  =  Q'(<r/T);  (3.9) 

this  follows  by  using  the  second  of  Eqs.  (3.2)  to  replace  Z  in  Eq.  (3.7).  This  relation  is 
important  because  Q'  is  positive  at  A  and  C  and  negative  at  B. 

The  character  of  the  eigenvalues  can  be  understood  using  these  foimulae  together  with 
Fig.  4.  In  addition  to  the  hyperbola  on  which  <7  =  0  and  the  parabola  on  which  Z  =  0 
[see  Eqs.  (3.2)],  Fig.  4  shows  the  hyperbola  on  which  Discrm  J  vanishes  [see  Eq.  (3.8)].  We 
draw  the  following  conclusions: 

(1)  Tr  J  <  0  at  all  critical  points; 

(2)  Det  J  >  0  at  A  and  C,  while  Det  J  <  0  at  B;  and 

(3)  Discrm  J  >  0  at  .4  and  B ,  whereas  Discrm  J  can  be  of  either  sign  at  C.  (For  “most” 
values  of  e  and  T,  Discrm  J  <  0  at  C;  in  particular,  Discrm  J  <  0  if  C  is  the  only 
critical  point.  But  it  is  possible  for  Discrm  J  to  be  positive  if  T  is  sufficiently  close  to 

Tm.) 

Standard  theory  of  nonlinear  planar  dynamical  systems  (see,  e.g.,  Ref.  [3,  Chap.  15])  now 
establishes  the  local  characters  of  the  critical  points  A ,  B,  and  C  in  Prop.  3.1: 

Proposition  3.2: 

(1)  A  is  an  attracting  node  (called  the  classical  attractor); 

(2)  B  is  a  saddle  point; 

(3)  C  is  either  an  attracting  spiral  point  or  an  attracting  node  (called  the  spurt  attractor). 
The  next  task  is  to  determine  the  global  structure  of  the  orbits  of  system  (3.1).  In 

this  direction,  we  employ  an  argument  based  on  a  suggestion  by  A.  Coppel  [4]. 


(3.7) 


(3.8) 


/-j(Z  + 1+e) 

V  ?(2#-1)  -1  ) 
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Proposition  3.3:  System  (3.1)  has  neither  periodic  orbits  nor  separatrix  cycles. 

For  the  proof,  notice  that  b  =  —  cr  ^  0  along  the  line  o  =  T,  so  that  no  periodic  orbit 
can  cross  this  line.  Since  all  critical  points  lie  in  the  half-plane  cr  <  T,  any  periodic  orbit 
(which  must  surround  some  critical  point,  by  the  Poincare-Bendixson  theorem)  necessarily 
lies  in  this  half  plane.  Let  the  scaled  time  variable  s  be  defined  by  ds  =  (T  —  cr)dt;  then 
with  a  prime  denoting  differentiation  with  respect  to  s,  system  (3.1)  transforms  to 

y  =  e~\Z  +  1)  -  <x/(T  -  a)  =:  p(a,  Z)  , 

E'  =  -c"V  -  Z/( T  -  a)  =:  g(cr,  Z)  . 

Since 

dp/ da  +  dq/dZ  =  — (2T  -  a)/(T  -  a)2  , 

is  strictly  negative  for  cr  <  T,  existence  of  periodic  orbits  and  separatrix  cycles  is  excluded 
by  Bendixson’s  criterion  [1]. 

To  understand  the  global  qualitative  behavior  of  orbits,  we  construct  suitable  invariant 
sets.  In  this  regard,  a  crucial  tool  is  the  identity 

d  {a-  +  (Z  +  l)2}  =  -2  [<72  +  (Z  +  ±)2  -  i]  ,  (3.12) 

which  is  obtained  by  multiplying  the  first  of  Eqs.  (3.1)  by  a  and  adding  the  second, 
multiplied  by  Z  +  1.  Thus  the  function  V((T,Z )  :=  <r2  -f  (Z  4-  l)2  serves  as  a  Lyapunov 
function  for  the  dynamical  system.  Notice  that  identity  (3.12)  is  independent  of  T  and  e. 

Let  r  denote  the  circle  on  which  the  right  side  of  Eq.  (3.12)  vanishes.  Also  let  Cr 
denote  the  circle  of  radius  r  centered  at  cr  =  0  and  Z  =  —  1,  along  which  V(cr,  Z)  —  r;  thus 
each  Cr  is  a  level  set  of  V.  The  circles  T  and  C\  are  shown  in  Fig.  5,  which  corresponds 
to  the  case  of  a  single  critical  point,  namely  the  spiral  point  C;  and  in  Fig.  7,  which 
corresponds  to  the  case  of  three  critical  points.  For  later  convenience,  denote  by  D  the 
point  where  C\  intersects  the  parabola  on  which  Z  =  0.  Also  notice  that  Eq.  (3.12)  implies 
that  the  critical  points  of  system  (3.1)  lie  on  T. 

If  r  >  1,  T  lies  strictly  inside  Cr ■  Consequently,  Eq.  (3.12)  shows  that  the  dynamical 
svstem  (3.1)  flows  inward  at  points  along  Cr-  Thus  the  interior  of  Cr  is  a  positively 
invariant  set  for  each  r  >  1.  Furthermore,  the  closed  disk  bounded  by  Ci,  which  is  the 
intersection  of  these  sets,  is  also  positively  invariant.  This  establishes. 

Proposition  3.4:  Each  closed  disk  bounded  by  the  circle  Cr  with  r  1  is  a  positively 
invariant  set  for  the  system  (3.1). 

The  above  results  will  be  used  to  determine  the  global  structure  of  the  orbits  of 
svstem  (3.1)  and  to  analyze  the  stable  and  unstable  manifolds  of  the  saddle  point  at  B. 


(3.10) 

(3.11) 
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Fig.  5:  The  phase  plane  when  the  spurt  attractor  C  is  the  only- 
critical  point. 

A.  Let  us  first  consider  the  structure  of  the  flow  when  there  is  a  single  critical  point, 
located  at  C;  see  Fig.  5.  As  shown  in  Prop.  3.2,  the  point  C  must  be  an  attracting  spiral 
point.  Because  there  are  no  periodic  orbits  (Prop.  3.3),  Prop.  3.4  implies  that  the  orbit 
through  each  point  in  the  phase  plane  tends  toward  C  as  t  approaches  infinity. 

In  the  application  to  the  shear  flow  problem,  we  are  interested  in  the  particular  solution 
of  Eq.  (3.1)  with  initial  data  <r  =  0  and  Z  =  0  (point  0 ).  The  orbit  of  this  solution  initially 
remains  inside  the  region  7Z  bounded  by  ODVCO  in  Fig.  5,  eventually  exits  through  the 
arc  CV  of  the  parabola,  and  finally  spirals  toward  C.  Indeed,  the  flow  is  directed  into 
71  along  OD  (where  V  <  0),  DV  (where  Z  =  0  and  a  <  0),  and  CO  (where  V  =  0  and 
&  >  0).  Therefore  the  orbit  through  0  must  leave  7Z  along  the  arc  VC  of  the  parabola, 
whereupon  it  spirals  into  C .  This  solution  is  illustrated  in  Fig.  6;  note  that  this  solution 
does  not  tend  to  a  periodic  orbit,  according  to  Prop.  3.3. 
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Fig.  6:  The  orbit  through  origin  when  the  spurt  attractor  C  is  the 
only  critical  point. 


B.  Next,  consider  the  case  when  there  are  three  critical  points,  illustrated  in  Fig.  7.  By 
Prop.  3.3.  there  are  no  periodic  orbits  or  closed  separatrix  cycles.  Propositions  3.2  and  3.4 
imply  that  the  orbit  through  any  point  in  the  plane  either:  tends  to  A;  tends  to  C ,  or 
tends  to  B  along  its  stable  manifold  as  t  approaches  infinity. 

For  a  more  detailed  analysis,  we  first  prove  that  the  closed  set  T  bounded  by  the 
curved  triangle  ODA  in  Fig.  7  is  positively  invariant  with  respect  to  Eqs.  (3.1).  This  is 
because  the  flow  is  directed  into  I  along  OD  (where  V  <  0),  DA  (where  Z  =  0  and  &  <  0), 
and  AO  (where  V  =  0  and  a  >  0).  One  consequence  of  the  invariance  of  I  is  that  the 
orbit  of  the  solution  of  Eq.  (3.1)  with  initial  data  cr  =  0  and  Z  =  0  flows  into  the  classical 
attractor  A  as  t  approaches  infinity. 

The  next  task  is  to  study  the  stable  manifold  for  the  saddle  point  B  in  Fig.  7.  Through 
this  point  we  have  drawn  the  circle  Cq  centered  at  <J  =  0  and  Z  —  —  1,  which  intersects  the 
parabola  at  E .  Let  <5  denote  the  closed  set  bounded  by  the  curved  triangle  ECB  in  Fig.  7. 
At  points  on  its  boundary,  the  flow  is  directed  as  follows  for  increasing  t:  outward  from  5 
along  EB  because  V  >  0;  outward  along  CB,  where  V  =  0  and  &  >  0;  and  inward  along 
CE  because  Z  =  0  and  a  >  0.  As  a  result,  one  branch  of  the  stable  manifold  at  B  must 
enter  5  through  the  arc  EC  and  remain  in  S.  Furthermore,  this  branch  must  have  crossed 
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Fig.  7:  Invariant  regions  in  the  case  of  three  critical  points. 

the  axe  CB'  of  T  (where  B'  is  the  reflection  of  B  about  the  Z- axis),  along  which  V  =  0 
and  a  ■  &  <  0,  because  it  could  not  have  crossed  the  axe  B' E  of  Cg.  The  other  branch  of 
the  stable  manifold  enters  B  through  the  sector  exterior  to  T  and  the  circle  Cg.  Having 
come  from  outside  of  the  circle  T,  where  V  <  0,  Prop.  3.4  implies  that  both  branches  of 
the  stable  manifold  must  have  come  from  infinity  in  the  phase  plane. 

For  the  purpose  of  analyzing  the  spurt  phenomenon  in  Sec.  5,  we  summarize  these 
results  as  follows. 

Proposition  3.5:  The  bzisin  of  attraction  of  A,  i.e.,  the  set  of  points  that  Sow  toward 
A  as  t  — *  oo,  comprises  those  points  on  the  same  side  of  the  stable  manifold  of  B  as  is 
A;  points  on  the  other  side  axe  in  the  basin  of  attraction  of  C .  Moreover,  the  arc  of  the 
circle  T  through  the  origin,  between  B  and  its  reSection  B'  is  contained  in  the  basin  of 
attraction  of  A. 

To  demonstrate  the  last  statement,  observe  that  the  flow  is  directed  into  the  region  bounded 
by  the  the  following  curves  in  Fig.  7:  C\  between  0  and  D\  the  parabola  between  D  and 
A;  T  between  A  and  B ;  Cg  between  B  and  B and  T  between  B'  and  O.  Therefore  this 
region  is  positively  invariant.  In  particular,  the  stable  manifold  for  B  cannot  cross  its 
boundary,  so  that  it  cannot  cross  T  between  B  and  B' . 

Finally,  consider  the  unstable  manifold  of  the  saddle  point  B.  Let  U\  be  the  set 


-18- 


bounded  by  the  axes  of  the  parabola  Z  =  0  and  the  hyperbola  a  =  0  in  Fig.  7,  between 
the  critical  points  B  and  A.  Along  the  open  arc  of  the  parabola  BA,  a  >  0,  while  along 
the  open  arc  of  the  hyperbola,  Z  >  0.  Therefore,  one  branch  of  the  unstable  manifold  at 
B  lies  in  Ui  and  connects  B  to  A.  Now  consider  the  set  Ui  bounded  by  BVC,  as  shown 
in  Fig.  7.  The  flow  is  directed  into  Ui  along  BV  and  CB,  so  the  second  branch  of  the 
unstable  manifold  at  B  remains  in  Ui  until  it  exits  through  the  arc  VC.  If  C  is  a  spiral 
point,  this  branch  enters  and  leaves  U 2  infinitely  often  as  it  spirals  into  C,  while  if  C  is  an 
attracting  node,  it  does  not  reenter  as  it  tends  to  C. 

To  summarize  the  above  description  of  the  dynamics  of  the  system  (3.1)  with  three 
critical  points,  the  reader  is  referred  to  Fig.  8,  which  shows  the  case  when  C  is  a  spiral 
point. 


Fig.  8:  Phase  portrait  in  the  case  of  three  critical  points,  with  C 
being  a  spiral. 


C.  Last  is  the  case  of  a  single  critical  point  at  A,  which  according  to  Prop.  3.2  is  an 
attracting  node.  Again,  because  there  are  no  periodic  orbits,  Prop.  3.4  implies  that  all 
orbits  are  attracted  to  A.  In  particular,  the  orbit  through  the  origin  remains  in  a  region 
that  is  analogous  to  the  region  X  in  Fig.  7. 
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4.  Phase  Plane  Analysis  for  a  Model  with  Two  Relaxation  Times 

The  purpose  of  this  section  is  to  study  the  qualitative  features  of  solutions  of  the 
system  ( J  SO  2 )  when  a  =  0.  We  focus  on  a  special  choice  for  the  parameters  defining 
the  model:  0  <  a\  <  1  and  a2  —  1  (so  that  5  =  0).  The  physical  significance  of  this 
choice  is  discussed  in  Sec.  6.  If  we  identify  k  with  e,  then  the  steady  states  for  this  model 
coincide  with  those  of  system  ( JSO )  [see  Eqs.  (2.32)-(2.34)];  again,  we  assume  £  <  1/8. 
Thus  the  second  relaxation  time  in  system  ( JS02 )  replaces  the  Newtonian  viscosity  in 
system  (JSO),  and  techniques  introduced  in  Sec.  3  can  be  applied.  The  question  we 
address  is:  to  what  extent  do  the  dynamics  of  the  present  model  resemble  the  dynamics 
of  system  (3.1)? 

When  a  =  0  and  f  —  f,  the  momentum  equation  in  system  ( JSO2 )  can  be  integrated 
to  show  that  the  total  shear  stress  T  :=  <7\+cr 2  equals  T(x)  —  —  fx;  equivalently,  cr1+&2  = 
0.  As  a  result,  adding  the  stress  equations  in  system  (JSO 2)  yields 

(Z\  4- 1  -j-  Z2  +  ,5)Ui  =  (?\  +  £  lj3(T  —  ai).  '  (4.1) 


Using  this  equation  to  eliminate  vx  leads  to  an  autonomous  system  of  three  ordinary 
differential  equations: 


(7  ,  n  T-(l 

crl  =  (Zl  +  l)7— 


e/3  l[Z  1  +  1  +  Z2  +  P\ 

- Z 1  , 


-  <7 1 


z  =-  T-(l -£/?->  1 
'  aie0-'[Zl  +  I  +  Z2+/)) 

*  <,!T  1  1  r-iR7 

^ ~  l)£/J-‘[Zi  +  l  +  Z2+/?|  0Z2' 


(4.2) 


The  assumption  5  =  0,  coupled  with  the  requirement  that  Z2  be  finite  as  t  — ►  —  cc, 
shows  that  Z2  =  0;  thus  the  last  of  Eqs.  (4.2)  can  be  eliminated.  Since  no  confusion  will 
arise  by  omitting  subscripts,  we  write  the  dynamical  system  as 


Notice  that  the  identity 

d_ 

dt 


(Z  +  1h/3-‘(Z  +  1  +  /J] 

T  -  (1  -  e0~l)<r 
Z  e(3-‘[2  +  l  +  ^l  Z' 


(ff2  +  (z  +  i)2}  =  -2y  +  (z  +  |)J-i] 


(4.3) 


(4.4) 


holds  for  system  (4.3)  just  as  it  does  for  system  (3.1). 

The  relationship  to  system  (3.1)  is  highlighted  by  scaling  the  time  parameter:  define 
the  scaled  time  s  to  satisfy  dt  =  ,3~l[Z  +  1  +  0]  ds,  and  let  a  prime  denote  the  derivative 
d/ds.  Then 


<r'  =  (z  +  1)(^t£) 

Z‘  =  -c  (^lS)-Z-0-'  {^+(Z+i)2-l} 


(4.5) 
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and 


d{^+(Z  +  l)2}  =  -2r1(Z  +  l  +  /31[<T2  +  (Z  +  |)2-i]  .  (4.6) 

Another  useful  way  of  writing  svstem  (4.3)  is  obtained  by  multiplying;  Eqs.  (4.3)  bv 
0~l[Z  +  1  +  0]  and  using  Eq.  (4.4): 

0  l[Z  +  1  +  p\&  =  ( Z  4-  1) 

-0~xab  +  Z  =  -a  ~  2  • 

It  follows  from  this  formulation  that  the  critical  points  satisfy  the  same  algebraic  system 
as  does  system  (3.1),  namely 


( Z  +  1  4-  s) 


+  £  —  0  , 


—  =  (= 
e  T\T 


-2  =  0. 


(4.3) 


Just  as  in  Prop.  3.1,  there  are  up  to  three  critical  points,  labelled  A,  B,  and  C,  for  each 
position  x  in  the  channel. 

The  character  of  the  critical  points  is  determined  by  the  Jacobian  J  of  the  right  side  of 
the  system  (4.3),  evaluated  at  the  critical  point.  Equivalently,  this  matrix  can  be  calculated 
using  Eqs.  (4.5)  or  Eqs.  (4.7),  which  imply  that 


0~l[Z+l+0]J 

/  -1(Z  +  l  +  e)  -?(f-l)  \ 

\  ?  (2f  -  l)  -  20- ‘a  -l-0-1(2Z  +  l)J 

_/  1  0  \  (- 7(Z  +  l  +  e) 

0-'{Z  +  l+0\J  ^ 


(4.9) 

(4.10) 


The  character  of  the  eigenvalues  of  J  can  be  determined  from  the  signs  of  the  trace 
of  J,  which,  according  to  Eq.  (4.9),  is  given  by 

—sd~l[Z  -f  1  +  /?]  Tr  J  =  Z  +  1  +2s  +  £0~1(2Z  4-  1)  ;  (4.11) 

the  determinant  of  J,  given  by 

— *•» 

£/3-1[Z  +  l+/3!DetJ  =  Z  +  l+«  +  ^-(^2j-l)  (=  -  l)  .  (4-12) 
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as  follows  from  Eq.  (4.10);  and  the  discriminant  of  J,  given  by 
er  j3~2  [Z  +  1  +  (3 ]2  Discrm  J 

=  [z  +  i-^-i(2z  +  i)]2-4T2(=-i)(2[i-£/3-‘]|-i)  .  (4'13) 

Because  Z  >  —  1  at  a  critical  point,  Tr  J  <  —2  +  f3~lm,  thus  Tr  J  <  0  at  all  critical 
points  provided  that  we  assume  that  (3  >  4.  Furthermore,  Eq.  (4.12)  shows  that  Det  J 
has  the  same  sign  at  the  critical  points  as  for  system  (3.1),  i.e.,  positive  at  A  and  C  and 
negative  at  B.  To  determine  the  signs  of  Discrm  J  at  the  critical  points,  first  notice  that 
it  is  positive  at  B  because  the  determinant  is  negative.  Furthermore,  Eq.  (4.13)  defines 
a  hyperbola  similar  to  that  defined  by  Eq.  (3.8)  and  drawn  in  Fig.  4;  it  crosses  the  line 
Z  =  —(1  —  £(5-1)/(  1  —  2s{3~1)  at  a  —  4(1  —  sf3~1)~1T  and  cr  =  T.  Assuming  that 
s3~l  <  4,  it  follows  that  the  discriminant  is  positive  also  at  A  because  a  b  <  a  a  <  T  and 
Zb  <  Z.\  <  0.  Similarly,  the  sign  of  the  discriminant  at  C  is  negative  unless  T  is  near  Tm, 
just  for  system  (3.1). 

As  in  Prop.  3.2,  these  results  determine  the  local  character  of  orbits  of  the  quadratic 
system  (4.3)  [and  the  equivalent  system  (4.5)]  near  critical  points:  A  is  an  attracting  node; 
B  is  a  saddle  point;  and  C  is  either  an  attracting  spiral  point  or  an  attracting  node.  We 
remark  that  a  spurt  attractor  C  that  is  a  node  for  system  (3.1)  might  be  a  spiral  point  for 
system  (4.3).  Furthermore,  Eq.  (4.4)  and  the  discussion  leading  to  Prop.  3.4  shows  that 
the  disk  bounded  by  the  circle  Ci  of  radius  1  centered  at  cr  =  0  and  Z  —  —  1  is  positively 
invariant  for  system  (4.3)  (cf.  Fig.  7).  We  assume  that  (3  >  1  so  that  the  line  Z  -f  l+f3  =  0 
does  not  cross  C\. 

To  rule  out  the  existence  of  periodic  orbits  and  separatrix  cycles  for  system  (4.5),  the 
analogue  of  Prop.  3.3,  consider  the  set  of  points  Q  inside  C\  such  that  cr  <  T(1  —  e/3~l)~1. 
Because  this  sec  contains  all  critical  points,  and  because  it  is  positively  invariant  [because 
d  =  —cr  0  along  the  line  cr  =  T(1  —  c/5-1)-1],  any  periodic  orbit  or  separatrix  cycle  must 
stay  within  Q.  Now  introduce  the  scaled  time  variable  u,  defined  by  du  =  (T  —  cr)  ds ,  into 
Eq.  (4.5)  to  obtain  the  equations  da / du  =  p(a.  Z)  and  dZ/du  =  q(a ,  Z).  A  straightforward 
calculation  shows  that 

dp/ da  +  dq/dZ  =  — (2T  -  a)/(T  -  a)2  -  (3~\2Z  +  1)/(T  -  a)  .  (4.14) 

This  quantity  is  strictly  negative  because  a  <  2 T  and  Z  >  —2  along  an  orbit  in  region  Q. 
(assuming  that  s3~l  <  ^).  Consequently,  Bendixson’s  criterion  excludes  periodic  orbits 
and  separatrix  cycles  for  the  system  (4.5)  and  hence  system  (4.3). 

The  remaining  analysis  of  the  phase  plane  for  the  system  (4.3)  is  almost  identical  to 
that  for  the  system  (3.1).  The  analysis  of  Sec.  3  relied  on  calculating  the  direction  of  flow 
across  five  curves:  the  circles  Cx  and  Ca,  along  which  V  is  constant  and  the  sign  of  V 
was  determined:  the  circle  V  where  V  =  0  and  a  had  a  definite  sign:  the  v  yperbola  where 
d=0  and  the  sign  of  Z  was  known;  and  the  parabola  where  Z  =  0.  For  system  (4.3),  V 
is  governed  by  the  same  equation  as  in  system  (3.1).  Furthermore,  Eqs.  (4.7)  shows  that 
(7  =  0  along  the  same  hyperbola,  and  that  the  sign  of  Z  on  this  hyperbola  is  the  same  as 


before.  Therefore  the  only  difference  concerns  the  parabola  where  =  e~l<j(T  —  a)  —  Z 
vanishes  and  Z  =  However,  =  [(2  +  e&~l)cr  -  T]  a  when  $  =  0,  so  that 

for  a  >  i T  (i.e.,  to  the  right  of  the  vertex  V),  $  has  the  same  sign  as  when  f3~l  =  0. 
Therefore  the  parabola  works  as  a  boundary  just  as  in  Sec.  3. 

Summary:  All  qualitative  features  of  the  dynamics  of  system  (3.1)  (except  possibly 
whether  C  is  a  node  or  a  focus)  carry  over  to  system  (4.3). 


5.  Qualitative  Features  of  ( JSO )  Based  on  Phase  Plane  Analysis 

The  principal  motivation  of  the  analysis  that  follows  is  to  explain  recent  numerical 
simulations  of  (JSO)  described  in  Refs.  [9,  11].  These  exhibited  several  effects  related 
to  spurt:  latency,  shape  memory,  and  hysteresis.  Fig.  9  shows  the  result  of  simulating  a 
“quasi-static”  loading  sequence  in  which  the  pressure  gradient  /  is  increased  in  small  steps, 
allowing  sufficient  time  between  steps  to  achieve  steady  flow  [9].  The  loading  sequence  is 
followed  by  a  similar  quasi-static  unloading  sequence,  in  which  the  driving  pressure  gradient 
is  decreased  in  steps.  The  initial  step  used  zero  initial  data,  and  succeeding  steps  used  the 
results  of  the  previous  step  as  initial  data.  The  resulting  hysteresis  loop  includes  the  shape 
memory  predicted  by  Hunter  and  Slemrod  [7]  for  a  simpler  model  by  a  different  approach. 
The  width  of  the  hysteresis  loop  at  the  bottom  can  be  related  directly  to  the  molecular 
weight  of  the  sample  [9]. 

In  this  section,  we  explain  spurt,  shape  memory,  and  hysteresis  using  the  results  of 
the  phase  plane  analysis  of  the  dynamical  system  (3.1)  corresponding  to  a  single  relaxation 
time;  the  discussion  can  also  be  based  on  the  system  (4.3)  for  two  relaxation  times.  Latency, 
which  occurs  for  (JSO)  when  £  is  sufficiently  small,  will  be  explained  in  Sec.  6  by  means  of 
an  asymptotic  analysis.  The  latter  will  also  be  used  in  Sec.  6  to  show  that  systems  (JSO) 
and  (JSO2)  exhibit  similar  features  when  s  is  sufficiently  small. 

We  consider  experiments  of  the  following  type:  the  flow  is  initially  in  a  steady  state 
corresponding  to  a  forcing  /0,  and  the  forcing  is  suddenly  changed  to  /  =  /0  4-  A/.  We 
call  this  process  “loading”  (resp.  “unloading” )  if  A /  has  the  same  (resp.  opposite)  sign  as 
f0.  The  initial  flow  can  be  described  by  specifying,  for  each  channel  position  x,  whether 
the  flow  is  at  a  classical  attractor  A  (x  is  a  “classical  point”)  or  a  spurt  attractor  C  (x  is 
a  “spurt  point”)  for  the  system  (3.1)  with  T  =  — / 0x.  We  shall  say  that  any  point  lying 
on  the  same  side  of  the  stable  manifold  of  B  as  is  A  lies  on  the  “classical  side”;  points 
lying  on  the  other  side  are  said  to  be  on  the  “spurt  side.”  The  outcome  of  the  experiment 
depends  on  the  character  of  the  phase  portrait  with  T  =  —  fx.  To  determine  this  outcome, 
we  need  only  decide  when  a  classical  point  becomes  a  spurt  point  or  vice  versa. 

To  understand  how  critical  points  move,  consider  Fig.  10.  First  notice  that  the  critical 
points  for  Eq.  (3.1),  for  any  value  of  T,  lie  on  the  circle  I\  which  is  independent  of  T. 
Let  M  =  (cr.v/,Z.vf)  denote  the  degenerate  (double  root)  critical  point  that  occurs  when 
T  =  T xi ,  i.e.,  at  “top  jumping”  in  Fig.  2;  and  let  m  =  (crm,Zm)  denote  the  degenerate 
critical  point  for  T  =  Tm ,  i.e.,  for  “bottom  jumping.”  These  points,  together  with  their 
reflections  A/'  =  (-cr.v/,  Z,w)  and  m'  =  (-<7m,Zm),  serve  to  divide  T  into  arcs:  the 

upper  arc  of  T  between  M  and  A/';  Tcj  the  lower  arc  between  m  and  m'\  and  T s,  the 
remaining  two  arcs  where  Z  €  [Zm,Z,vf].  For  any  value  of  T,  positive  or  negative,  the 
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Fig.  9:  Hysteresis  under  cyclic  load:  normalized  throughput  6Q 
vs.  wall  shear  stress  Inwall  [9j- 

classical  attractor  A  lies  in  T.4,  the  spurt  attractor  C  lies  in  Tc,  and  the  saddle  point  B 
lies  in  Y g.  Furthermore,  as  \T\  is  increased,  the  critical  points  A  and  C  move  downward 
along  T,  while  B  moves  upward.  This  follows  from  Eq.  (3.3)  by  differentiating  the  relation 
Q(cr/T)  =  0,  which  determines  how  cr/T  varies  with  T,  and  by  using  the  first  of  Eqs.  (3.2). 

The  principle  mathematical  properties  of  the  dynamical  system  (3.1)  that  determine 
the  outcome  of  loading  and  unloading  experiments  are  embodied  in 

Proposition  5.1: 

(1)  A  classical  point  Ar>  for  the  initial  forcing  f0  lies  in  the  domain  of  attraction  of  the 
classical  attractor  .4  for  f,  provided  that  A  exists  (i.e.,  \fx\  <  T m); 

(2)  A  spurt  point  Cq  for  the  initial  forcing  f0  lies  in  the  domain  of  attraction  of  the  spurt 
attractor  C  for  f  unless  (a)  C  does  not  exist  (i.e.,  \fx\  <  Tm);  or  (b)  C  lies  on  the 
classical  side  of  the  stable  manifold  of  the  saddle  point  B  for  f. 

This  follows  because,  as  shown  in  Secs.  3  and  4  (see  Prop.  3.5),  the  stable  manifold  of 
B.  which  separates  the  domain  of  attraction  of  .4  from  that  of  C,  never  crosses  the  arc 
rv  We  remark  that  (2b)  occurs  only  when  |A/|  is  large  enough.  In  particular,  it  never 
occurs  during  a  quasi-static  processes,  in  which  the  flow  is  allowed  to  attain  equilibrium 
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Fig.  10:  Location  of  critical  points. 

after  each  of  many  small  increments  in  forcing. 

A.  Startup 

As  a  first  experiment,  consider  starting  with  /0  =  0  and  loading  to  /  >  0.  Thus 
the  initial  state  for  each  x  lies  at  the  origin  a  =  0,  Z  =  0.  Then  according  to  5.1(1) 
above,  each  x  €  [—1/2.0]  such  that  /|xj  <  is  a  classical  point,  while  the  x  for  which 
f\x\  >  T \i  are  spurt  points  (because  there  is  no  classical  attractor).  Consequently,  we 
draw  two  conclusions: 

Proposition  5.2: 

(a)  If  the  forcing  is  subcriticaJ  (i.e.,  f  <  /crit  :=  2 T m),  the  asymptotic  steady  Sow  is 
entirely  classical. 

(b)  If  the  forcing  is  supercritical  (f  >  /criJ,  there  is  a  single  kink  in  the  velocity  proSle 

( see  Fig.  3),  located  at  x,  =  those  x  6  [— 1/2,  x.),  near  the  wall ,  are  spurt 

points,  whereas  x  6  (x.,0],  near  the  centerline,  are  classical. 

The  solution  in  case  (b)  can  be  described  as  “top  jumping”  because  the  stress  T«  =  T \i 
at  the  kink  is  as  large  as  possible,  and  the  the  kink  is  located  as  close  as  possible  to  the 
wall. 


B.  Loading 

Next.  consider  increasing  the  load  from.  /0  >  0  to  /  >  /0.  A  point  x  that  is  classical 
for  f0  remains  classical  for  /  unless  there  is  no  classical  attractor  for  T  =  —fz,  he., 
f'x\  >  T A  spurt  point  x  for  f0.  on  the  other  hand,  is  always  a  spurt  point  for  /.  This 
is  true  because:  (a)  the  spurt  attractor  always  exists  at  the  larger  total  stress  T  =  —fx; 
and  (b '  the  stable  manifold  of  B  does  not  cross  the  arc  of  T  between  C  and  B ,  and  Cq  lies 
on  this  arc.  so  that  C0  lies  on  the  spurt  side  of  the  stable  manifold  of  B. 

As  a  result,  a  point  in  x  in  the  channel  can  change  only  from  a  classical  attractor  to 
a  spur:  attractor,  and  then  only  if  f\x\  exceeds  T When  /  is  chosen  to  be  supercritical, 
loading  causes  the  position  x .  of  the  kink  in  Fig.  3  to  move  away  from  the  wall,  but  only 
to  the  extent  that  it  must:  a  single  jump  in  strain  rate  occurs  at  x .  =  — T\i/f ,  where 
the  total  stress  is  T.  =  T \j.  These  conclusions  are  valid,  in  particular,  for  a  quasi-static 
process  of  gradually  increasing  the  load  from  fQ  =  0  to  /  >  /crit. 

C.  Unloading 

Now  consider  unloading  from  /0  >  0  to  /  <  fQ;  assume,  for  the  moment,  that  f  is 
positive.  Here,  the  initial  steady  solution  need  not  correspond  to  top  jumping.  For  this 
type  o:  unloading,  a  point  x  that  is  classical  for  /0  always  remains  classical  for  /:  the 
classical  attractor  for  /  exists  because  f\x\  <  fQ\x |.  By  contrast,  a  spurt  point  x  for  /0 
can  become  classical  at  /.  This  occursjf:  (a)  the  total  stress  T  =  —fx  falls  below  Tm;  or 
( b  >  the  spurt  attractor  Co  for  T  =  —f0x  lies  on  the  classical  side  of  the  stable  manifold 
of  the  saddle  point  B  for  T  =  —fx  [see  Prop.  5.1(2b)j.  The  latter  possibility  can  cause 
the  formation  of  a  second  kink  in  the  steady  velocity  profile,  which  arises  from  a  second 
jump  in  the  plot  of  steady  shear  stress  vs.  strain  rate  shown  in  Fig.  11.  A  new  region  of 
classical  flow  develops  next  to  the  wall,  and  the  internal  region  of  spurt  flow  lies  between 
the  two  classical  regions.  Of  course,  tills  occurs  only  if  | A/j  =  \f  —  f0 1  is  sufficiently  large, 
i.e.,  the  unloading  is  abrupt.  However,  it  is  impossible  to  have  classical  flow  regions  to 
start  spurting,  even  under  sudden  finite  loading,  unless  the  classical  attractor  disappears. 
In  a  quasi-static  experiment,  in  which  the  forcing  is  decreased  gradually,  spurt  points  can 
become  classical  only  if  /|xj  <  Tm. 

D.  Quasi-Static  Loading  and  Unloading:  Shape  Memory  and  Hysteresis 

We  now  discuss  the  quasi-static  loading-unloading  cycle,  illustrated  in  Fig.  9,  based 
on  results  of  Subsections  B.  and  C.  Specifically,  recall  that  there  axe  two  ways  in  which  a 
spurt  point  x  can  become  classical  during  unloading:  (i)  if  T*  decreases  to  Tm,  or  (ii)  if 
case  (bi  in  Subsection  G  occurs.  For  a  fixed  £,  the  latter  can  be  avoided  in  quasi-static 
unloading  by  keeping  the  increments  |A|  sufficiently  small. 

Consider  a  quasi-static  loading- unloading_cycle  that  starts  with  zero  forcing,  loads  to 
a  7max_>  7r.ru-  3X1(1  then  unloads  to  /final  <  /crit;  of  course,  /final  may  be  zero.  At  first, 
while  /  <  fcrn.  the  flow  is  classical  (from  the  origin  to  a  in  Fig.  9).  When  /  becomes 
supercritical  (from  a  to  6  in  Fig.  9),  a  kink  forms  at  the  wall;  the  total  stress  at  the  kink _i_s 
T.  =  T\f.  and  the  position  x,  moves  out  from  the  wall  according  to  the  rule  x,  =_— T^vr//. 
The  thickest  spurt  layer  occurs  at  /max  (b  in  Fig.  9),  for  which  x„  —  xmax  =  ~T \f / f  max. 
Upon  unloading  (from  b  to  c  to  d  in  Fig.  9),  the  layer  position  x  .^remains  fixed  at  xraax, 
and  the  stress  is  T.  =  (7/7ma «)7i/,  so  long  as  T.  remains  above  Tm;  this  occurs  between 
b  and  c  in  Fig.  9  and  corresponds  to  shape  memory,  in  the  sense  that  the  layer  position 
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Fig.  11:  Generation  of  flow  with  two  kinks,  giving  rise  to  an  inter¬ 
nal  spurt  layer.  The  sequence  of  relevant  strain  rates,  vx,  labeled 
1  -  5,  is  shown  as  the  stress,  w ,  varies  linearly  across  the  channel. 

remains  fixed,  and  points  on  the  interval  [— 1/2,  x«)  remain  spurt  points,  while  points  on 
the  interval  (x„,  0]  remain  classical.  Notice  that  this  behavior  gives  rise  to  hysteresis,  in  the 
sense  that  no  part  of  the  loading  curve  in  Fig.  9  is  retraced  until  the  entire  flow  becomes 
subcritical.  If  /final  is  small  enough,  then  eventually_T,  falls  below  Tm,  whereupon  the 
layer  position  x„  moves  toward  the  wall:  x«  =  — Tm//,  and  the  stress  at  the  kink  is  Xm; 
this  occurs  between  c  and  d  in  Fig.  9  and  corresponds  to  loss  of  shape  memory.  If  x, 
reaches  the  wall  (at  d  in  Fig.  9),  then  the  flow  becomes  entirely  classical. 

Summary:  A  kink  moves  away  from  the  wa 11  under  loading  if  it  is  a  top  jumper;  it  moves 
toward  the  wadi  under  unloading  if  it  is  a  bottom  jumper;  and  otherwise  it  remains  fixed. 
The  hysteresis  loop  opens  from  the  point  at  which  unloading  commences;  no  part  of  the 
unloading  path  retraces  the  loading  path  until  point  d  of  Fig.  9  is  reached. 
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E.  Reloading  and  Flow  Reversal 

It  is  possible,  of  course,  to  apply  more  complex  loading  and  unloading  sequences  than 
those  described  in  Subsection  D.  For  example,  a  loading  sequence  can  be  followed  by 
unloading,  which,  in  turn,  is  followed  by  reloading.  Arguments  extending  the  ones  given 
above  can  be  used  to  make  qualitative  predictions  in  such  cases.  In  particular,  Prop.  5.1 
implies  that  the  layer  position  can  remain  fixed  upon  reloading  (i.e.,  shape  memory  is 
retained). 

Next  we  consider  flow  reversals,  which  occur  when  fQ  and  /  have  opposite  signs. 
In  flow  reversal  with  loading,  i.e.,  |/|  >  |/0|,  a  mild  generalization  of  the  argument  of 
Subsection  B.  that  the  kink  position  x,  moves  away  from  the  wall,  and  all  conclusions 
hold  with  appropriate  change  in  sign.  To  see  this,  it  is  helpful  to  examine  Fig.  10  with 
appropriate  parabolas  and  circles  T  and  Cb  constructed  also  in  the  left  half-plane;  this 
will  also  be  useful  in  what  follows.  A  striking  phenomenon  has  been  observed  in  numerical 
experiments  using  data  in  Ref.  [17],  in  the  case  flow  reversal  in  which  a  spurted  solution 
is  unloaded,  i.e.,  f0  and  f  have  opposite  signs  and  |/|  <  |/0|.  Surprisingly,  numerical 
experiments  indicate  that  in  such  a  process  the  layer  position  can  remain  unchanged,  and 
shape  memory  can  be  retained.  For  this  to  happen,  the  magnitude  of  the  stress  at  the  layer 
boundary  must  not  fall  below  Tm.  Given  this  fact,  there  is  a  variety  of  phase  portraits 
constructed  with  the  aid  of  Fig.  10  that  can  lead  to  retention  of  shape  memory  in  flow 
reversal  with  unloading.  Perhaps  the  simplest  way  this  can  occur  is  to  load  continuously 
to  a  supercritical  load  level  fQ  >  0  to  obtain  a  top-jumping  solution  of  the  system  (3.1) 
in  which  all  spurt  points  x  €  [— 1/2,  x.)  have  only  a  single  critical  point  C0.  Reversing 
this  flow  with  unloading,  leads  to  solutions  in  which  all  points  x  £  (x,,0]  remain  classical; 
observe  that  circles  Cb  in  Fig.  10  are  symmetric  about  the  Z-axis,  and  by  symmetry, 
Prop.  5.1(1)  implies  the  result.  It  is  clear  that  points  x  £  [— 1/2,  x„)  for  which  T  >  T m 
will  remain  spurt  points.  The  only  remaining  issue  concerns  points  x  £  [— 1/2, x.)  for 
which  T m  <  T  <  T \f.  For  these  to  remain  spurt  points,  the  stable  manifold  of  B  must 
extend  fax  enough  into  the  right  half-plane  to  enclose  the  spurt  attractors  Co.  General 
criteria  for  this  situation  to  occur  are  not  so  easy  to  establish;  however,  it  will  be  seen  in 
Sec.  6  that  for  ;  as  small  as  in  Ref.  [17],  this  situation  is  likely  to  occur.  By  contrast,  there 
are  many  ways  to  lose  shape  memory  particularly  if  e  is  not  small. 


6.  Asymptotic  Analysis  of  ( JSO )  and  (JSO2)  for  a  =  0  and  Small  e 

The  analysis  of  system  (JSO)  leads  to  analytical  description  that  accurately  reflects 
the  experimental  results  of  Vinogradov  et  ai.  [17].  It  raises  a  question  of  the  physical  in¬ 
terpretation  of  the  Newtonian  viscosity.  This  was  chosen  to  give  the  total  stress  vs.  strain 
rate  curve  local  extrema  and  evidently  led  to  an  empirical  model  with  mathematical  fea¬ 
tures  that  reproduce  observed  behavior.  To  make  the  predictions  of  this  model  credible, 
the  empirical  parameter  s  must  be  associated  with  a  measurable  physical  quantity.  Since 
there  is  no  solvent  associated  with  the  undiluted  systems  studied  by  Vinogradov  et  a/., 
the  most  likely  interpretation  of  the  Newtonian  contribution  is  that  is  somehow  represents 
the  effects  of  shorter  relaxation  times.  At  first  sight  this  seems  puzzling  in  that  the  full 
system  (2.9)  is  not  classifiable  according  to  type  when  rj  ^  0  but  is  hyperbolic  (away  from 
the  elliptic  region)  when  77  =  0;  a  small  amount  of  Newtonian  viscosity  seems  to  make  a 
large  difference  in  the  mathematical  character  of  system  (2.9).  In  Sec.  4,  we  saw  that  in  the 
region  of  parameter  space  characteristic  of  the  data  of  Vinogradov  et  ai.,  system  (JSO)  and 
system  (JSO2)  exhibit  the  same  steady  behavior  and  very  similar  dynamic  behavior.  This 
suggests  that  it  may  be  valid  to  approximate  system  (JSO2)  by  the  simpler  system  (JSO). 
In  this  section,  we  show  that  the  dynamic  behavior  of  subcritical  solutions  and  of  critical 
solutions  during  the  latency  phase  in  the  two  models  agrees  to  0(e2).  Since  the  existence, 
location,  and  nature  of  the  critical  points  of  system  (JSO)  and  system  (JSO2)  is  identical 
(for  3  >  1)  and  the  phase  portraits  are  qualitatively  the  same,  the  only  possible  difference 
between  the  two  models  can  be  seen  to  lie  in  the  precise  timing  of  the  “Newtonian  phase” 
[11]  and  the  details  of  the  spurt  transition  following  latency. 

The  following  calculations  involve  system  ( JSO2 )  but  are  valid  for  system  (JSO)  in 
the  limit  as  (3  — >  00.  A  major  point  is  that  many  of  the  results  are  independent  of  (3.  We 
consider  the  case  when  e  •<  1;  this  is  true  for  the  data  of  Vinogradov  et  ai.that  exhibit 
spurt,  where  s  is  in  the  range  O(10~"2)  to  0(1O'“3).  Equation  (4.4)  shows  that  solutions 
of  system  (JSO 2),  with  initial  data  at  the  origin,  on  time  scales  of  order  e,  very  nearly 
follow  the  arcs  of  circles  defined  by  the  right  side  of  the  equation.  Equation  (4.5)  shows 
that  such  solutions  will  change  on  an  0(e)  time  scale  until  a  =  T  —  0(e).  The  time  during 
which  this  occurs  is  referred  to  as  the  “Newtonian  phase.”  In  principle,  it  is  possible  to 
obtain  higher  order  estimates,  but  for  our  purposes,  it  will  be  sufficient  to  omit  terms  of 
order  £  and  higher  and  approximate  the  value  of  Z  at  the  end  of  the  Newtonian  phase  by 
that  value  on  the  circle  corresponding  to  <7  =  T.  This  value  is 


There  are  two  important  points  to  be  made  here.  First,  for  T  >  1  spurt  starts  immediately 
and  there  is  no  Newtonian  or  latent  phase.  Second,  the  behavior  of  systems  (JSO)  and 
(J SOn)  differs  only  by  a  factor  of  3~l(Z  4-  1  +  j3)  in  the  time  scale  and,  to  0(e),  they  have 
the  same  orbits  during  the  Newtonian  phase.  We  now  consider  system  (JSO2)  during  the 
latent  phase,  which  begins  when  Eq.  (6.1)  holds  approximately. 

We  make  a  formal  asymptotic  expansion  of  the  solution  (cr,  Z)  of  the  system  (4.6)  in 
powers  of  e, 

<7  =  <?q  -4*  e<?i  +  e2 &2  +  . . . 

Z  =  Zq  4-  eZ\  +  e2Z2  +  . . . 


(6.2) 


Multiplication  of  Eqs.  (4.7)  by  e  and  substitution  of  Eqs.  (6.2)  readily  yields  the  following 
result  upon  comparing  terms  through  first  order  in  e: 


Co  =  T 
-  ? 

Zq  +  1  (6-3) 


Z .v  provides  an  initial  value  for  Zq  in  the  solution  of  the  third  equation  for  Zq.  Thus  cr 
stays  nearly  constant  while  Z  varies  at  zero  order.  This  is  a  key  point:  The  long  relaxation 
mode  response  of  the  system  is  a  rapid  development  of  a.  At  early  time,  during  latency,  a 
balances  the  pressure  gradient  by  assuming  the  linear  profile  of  T.  During  this  time,  the 
normal  stresses,  which  are  self-equilibrating  in  that  they  play  no  direct  role  in  balancing  the 
pressure  gradient,  grow  on  a  very  different  time  scale.  The  solution  of  the  third  equation  of 
(6.3),  which  we  call  the  “latency  equation”  is  independent  of  all  parameters  of  the  model 
except  T.  For  1/2  <  T  <1,  latency  ends  when  Zq  is  sufficiently  close  to  —1;  for  this  value 
of  Zq,  becomes  singular,  and  spurt  ensues.  For  T  <  1/2,  as  we  have  seen  from  the  phase 
plane  analysis,  orbits  starting  at  the  origin  do  not  move  to  the  spurt  attractor.  The  latency 
time  can  be  easily  found  by  solving  for  the  inverse  function  in  the  third  equation  of  (6.3), 
which  can  be  done  analytically;  the  duration  of  the  Newtonian  phase  can  be  neglected, 
and  the  integration  can  be  performed  from  Zq  =  Zn  to  Zq  —  —  1,  as  follows 


^latency 


/, 


-  —  — dZ  =  \  In 
zN  Z^  +  Z  +  T  2 


Zff  +  Zn  +  T 


1  . 
-1 - -=  =  -  tan  1 

T-  -  1 


tan  1 


\  Y  4T“  —  1/ 


(6.4) 


The  resulting  latency  time  is  non-dimensional  and  thus  scales  with  A-1  for  a  fixed  T.  The 
non-dimensional  latency  time  decreases  with  increasing  T,  because  the  initial  value  given 
by  Eq.  (6.3)  gets  closer  to  one  and  the  integrand  of  Eq.  (6.4)  decreases.  For  T  <  1/2, 
one  finds  that  the  integral  with  the  given  limits  does  not  exist,  so  that  in  this  sense,  the 
latency  time  of  a  classical  solution  is  infinite.  Note  that  these  results  for  system  (<750 2) 
are  independent  of  0,  and  thus  apply  equally  to  system  ( JSO ). 

Numerical  initial  values  for  the  coefficients  of  Eqs.  (6.2)  at  the  first  two  orders  can 
be  determined  from  matched  asymptotic  expansions  in  the  spirit  of  Ref.  [10].  It  turns 
out  that  all  expansion  terms  through  order  one,  except  Z\,  of  an  inner  expansion  with 
time  stretched  by  tm  :=  t/e  have  finite  asymptotic  values  as  tm  — *  cc.  A  value  of  can 
be  obtained  by  taking  the  value  obtained  when  the  coefficients  of  zero  order  are  within 
0(s 2)  of  their  asymptotic  values,  and  the  coefficients  of  first  order  axe  within  0(e)  of  their 
asymptotic  values.  This  procedure  yields  initial  data  for  the  outer  expansion,  Eq.  (6.2) 
from  the  inner  expansion;  it  is  not  sensitive  to  the  precise  values  of  the  “O-constants” 
involved. 
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Continuing  the  outer  expansion  at  second  order,  we  find  that  a2  and  Zx  satisfy 

f3  1(Z  +  l)& i  +  cr i  =  ~(Zq  +  1)^2  ~  ( Z\  +  l)oq 
~  (3  *<Jq&i  -f-  Z\  —  O'0<J2  cf  —  Z\  ^  ^ 

Since  ax  is  known,  the  first  equation  reduces  to  an  algebraic  equation  for  cr2  in  terms  of 
known  quantities  and  Z\,  substitution  into  the  second  equation  gives  an  ordinary  differ¬ 
ential  equation,  with  known  coefficients,  for  Z\.  There  is  an  interesting  feature  of  such  a 
solution.  Multiplying  the  first  equation  by  ctq  and  the  second  by  Zq  +  1  and  adding  yields 

(Zq  +  1  )Z\  +  CTqCTx  =  (Zq  +  1  )&!  —  (Z0  -f-  1  )Z\  +  0’q(Zi  +  l)^  (6-6) 

So  Z\  is  determined  by  an  ODE  whose  coefficients  are  independent  of  /?,  and  if  initial  data 
is  given  for  Zx  that  is  independent  of  (3 ,  then 


~  f\(!3,T ,  Zq,Z\) 
Zl=f2(T,Zo,Z1) 


Thus  systems  ( JS02 )  and  ( JSO )  agree  to  two  orders  in  s  in  their  behavior  during  the 
latent  phase,  when  their  initial  data  agree  to  that  level  of  accuracy.  However,  the  value 
of  Zx  obtained  from  the  matching  procedure  described  above,  and  the  precise  value  f,  at 
which  the  matching  criteria  are  met  are  functions  of  /?.  This  leads  to  the  conclusion  that 
the  two  systems  (JS02)  and  (JSO)  agree  to  order  0(e)  in  the  Newtonian  phase;  during 
the  latency  phase,  the  difference  between  the  models  accumulates  further  only  at  order 
0(e2).  Moreover,  changing  the  values  of  the  O -constants  in  matching  does  not  alter  this 
conclusion. 

By  a  perturbation  analysis  of  the  roots  of  system  (2.22)  for  small  s,  one  finds  that  at 
the  spurt  attractor 


a  = 


e 

T 


Z  =  - 


1  +  rrv  +  0(e 2) 
T~ 

[1  +  0(£2)] 


(6.8) 


An  important  consequence  of  (6.3)  relates  to  retention  of  shape  memory  in  flow  rever¬ 
sal.  Note  that  the  distance  from  Co  to  C  is  of  order  s,  and  recall  from  Sec.  5E  that  shape 
memory  is  lost  only  if  the  stable  manifold  of  B  separates  Cq  and  C . 

Continuing  the  asymptotic  analysis  of  spurt  dynamics,  inspection  of  Eqs.  (3.1)  or  (4.3) 
indicates  that  in  the  early  stages  of  spurt,  when  Z  is  not  near  —1  and  <r  has  deviated  from 
T  after  latency,  the  derivatives  of  the  system  are  0(e-1).  A  stretched-time  series  like 
that  for  the  Newtonian  phase  would  be  appropriate;  analysis  of  the  Newtonian  phase  then 
suggests  that  the  two  systems  could  differ  in  their  spurt  dynamics  at  least  at  first  order, 
but  of  course,  they  would  end  up  in  the  same  steady  states.  The  conclusion  is  that  systems 
(JS02)  and  (JSO)  have  the  same  critical  points  and  qualitatively  similar  orbits.  In  the 
asymptotic  range  studied  here,  they  have  virtually  identical  latency  times. 

The  remaining  task  is  to  study  possible  differences  between  (JSO)  and  (JS02)  that 
may  occur  during  the  spurt  phase,  and  predict  dynamic  behavior  that  may  be  observable 
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in  experiments.  The  analysis  that  follows  will  determine  the  detailed  structure  of  orbits 
spiraling  into  the  spurt  attractor  C,  for  example  in  Fig.  7,  in  as  large  a  neighborhood  of  C 
as  permitted.  One  objective  is  to  determine  the  frequency  at  which  orbits  spiral  toward  C. 
We  use  the  scaled  system  (4.5)  that  approximates  JSO 2.  Define  X  :=  (cr,  Z)T  and  write 
X  =  A'o  +  X,  where  XQ  =  (<7o,  Z0)T  is  a  critical  point.  Then  the  rescaled  system  (4.5)  is 
of  the  form  X'  =  F(X);  by  Taylor’s  theorem,  the  latter  can  be  written  as 

X'  =  3X  +  \xtD*F(X0)X  ,  (6.9) 

where  D  denotes  differentiation  with  respect  to  A,  F(Xq)  =  0,  DF(Xq)  =  J  of  Eq.  (4.9), 
and  D2F  is  the  Hessian.  The  Taylor  series  stops  at  second  order  because  the  system  is 
quadratic.  We  convert  system  (6.9)  to  polar  coordinates  in  the  phase  plane,  with  origin  at 
A’o,  bv  defining  r  and  9  implicitly: 


a  =  <?q  +  r  cos  9 
Z  =  Zq  4-  r  sin  9  , 


and  we  obtain 


(lnr)'  = 


ad'  + 


9'  = 


a'Z  -  aZ' 


(6.10) 


(6.11) 

(6.12) 


We  shall  evaluate  these  derivatives  centered  at  the  spurt  attractor,  C,  when  it  is  a  spiral 
point.  The  ratio  of  these  two  derivatives,  dlnr/d.9,  when  the  numerator  is  negative,  is 
the  rate  of  motion  toward  the  attractor  compared  with  the  angular  frequency  of  spiraling 
around  the  attractor. 

Applying  Eqs.  (6.9)  and  6.10  leads  directly  to 


X'  =  JX  + 


\e~lcr2  —  /?_1r2  )  • 


(6.13) 


Note  that  it  is  much  easier  to  determine  the  quadratic  terms  from  direct  inspection  of 
system  (4.5),  rather  than  formal  evaluation  of  the  second  derivative  in  Eq.  (6.9).  It  is 
a  simple  matter  to  use  Eq.  (6.13)  to  derive  a  general  expression  for  dlnr/d^,  but  the 
expression  is  rather  complicated;  it  is  greatly  simplified  when  evaluated  at  C,  using  the 
asymptotic  formulas  Eqs.  (6.8)  to  evaluate  cr0,  Z0.  Rescaling  back  to  real  time,  the  result 
is 


dlnr 


d9  t 


£  1 

+  ±2L  _  g-i 

(  '  r^T  * 

\z  -  ■ 

l 

z  T 

r* 

) 

T -a 

-e  | 

(i(i -£)+/»-> 

cr  — 

vZ 
r 3  _ 

) 

+  0(£2) 


(6.14) 


There  are  several  restrictions  on  the  use  of  Eq.  (6.14)  contained  in  the  inequalities  that 
follow.  These  are  asymptotically  valid  to  order  0(e);  that  is  they  may  require  the  addition 


-32- 


of  a  term  of  magnitude  at  most  order  0(c)  to  make  them  true.  The  first  restriction  is 
that  the  change  in  the  time  scale  (see  Eq.  (4.5))  be  valid  (i.e.,  that  the  denominators  in 
Eq.  (4.3)  do  not  vanish).  This  is  assured  by  requiring 

(3~lZ>- 1.  (6.15) 

Also  r  and  T  must  be  restricted  so  that  the  solution  spirals  into  the  spurt  attractor,  C, 
e.  g.,  by  requiring  the  numerator  on  the  right  side  of  Eq.  (6.13)  to  be  negative,  while 
the  denominator  is  positive.  This  yields  the  following  asymptotic  inequalities  that  are 
sufficient  to  produce  spiraling: 


0  <  r  <  r0  <  min{/?  —  1,T} 

f  >  - I - 

~  2  —  2j3~1(r0  +  1)  ‘ 


(6.16) 


These  inequalities  may  be  applied  by  restricting  r0_ on  the  basis  of  (3  —  1  alone  and  applying 
the  second  inequality  in  Eq.  (6.16)  to  determine  T;  if  necessary,  r0  can  then  be  further 
restricted  to  satisfy  the  first  of  inequalities  (6.16).  Note  that  T  can  never  be  smaller  than 
1/2.  so  that  C  is  always  a  spiral  point.  Spiral  orbits  are  observed  for  T  smaller  than  1/2, 
but  restrictions  on  r  that  are  more  optimal  than  made  in  (6.16)  would  be  required  in  such 
cases.  We  emphasize  that  the  restrictions  imposed  on  r  are  much  less  severe  than  would 
be  implied  by  retaining  only  linear  terms  in  Eq.  (6.9). 

The  key  observation  for  determining  the  frequency  of  spiraling  toward  C  is  that  both 
linear  and  nonlinear  terms  of  order  c-1  cancel  in  the  radial  equation,  Eq.  (6.11),  and  thus 
in  the  numerator  of  Eq.  (6.14);  this  can  be  seen  from  Eqs.  (4.9)  and  (6.13)  (note  there 
are  both  linear  and  nonlinear  contributions  to  the  angular  frequency  of  order  e  l).  Thus 
the  orbics  are  elliptical  spirals  with  decay  rates  of  order  0(e)  in  the  radial  component;  as 
e  — *  0,  the  orbits  tend  to  circles.  These  conclusions  hold  as  far  from  C  as  permitted  by  the 
inequalities  of  Eq.  (6.16).  We  can  also  see  that  for  (3  >  1,  the  effect  of  0  is  to  modify  the 
decay  rate  of  the  radial  component  to  relative  order  0(1),  but  the  effect  is  only  at  relative 
order  0(c)  in  angular  frequency.  For  small  values  of  £,  the  decay  rate  is  sufficiently  small 
that  the  behavior  could  be  mistaken  for  undamped  oscillation  of  the  stresses. 
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7.  Physical  Implications 

The  results  of  Sec.  5  provide  an  explanation  of  the  spurt  phenomenon  in  the  context 
of  a  model  that  is  analyzed  in  sufficient  detail  to  identify  a  number  of  new  features  that 
should  experimentally  observable.  We  will  summarize  these  predictions  below.  First,  it  is 
important  to  note  that  there  is  a  conclusion  of  our  analysis  that  on  the  accuracy  of  these 
predictions:  When  momentum  transfer  is  not  important  in  the  dynamic  process  (that  is, 
when  it  is  a  valid  approximation  to  set  a  to  zero  in  ( JSO )  and  ( JS02 ),  a  small  amount  of 
Newtonian  viscosity  in  system  {JSO)  accurately  mimicks  the  behavior  of  system  ( JS02 ) 
with  5  =  0  and  small  £.  This  is  true  in  spite  of  the  fact  that  ( JSCN )  is  hyperbolic  while 
{JSO)  cannot  be  classified  according  to  type.  We  expect  that  this  observation  extends  to 
other  models  and  provides  some  quantitative  reinforcement  to  the  belief  among  rheologists 
that  Newtonian  viscosity  can  be  a  useful  modeling  tool  in  undiluted  systems. 

One  of  the  widely  accepted  explanations  of  spurt  and  similar  observations  is  that  the 
presence  of  the  wall  affects  the  dynamics  of  the  polymer  system  near  the  wall.  Conceivably, 
there  could  be  a  variety  of  “wall  effects,”  the  most  obvious  is  the  loss  of  chemical  bond 
between  wall  and  fluid,  or  wall  slip  [5].  Perhaps  the  most  distinguishing  feature  of  our 
alternative  approach  is:  it  predicts  that  spurt  steins  from  a  material  property  of  the 
polymer  and  is  not  related  to  any  external  interaction.  The  spurt  layer  forms  at  the  wall 
in  situations  such  as  top  jumping  because  the  stresses  are  higher  there;  for  the  same  reason, 
of  course,  is  chemical  bonds  would  break  at  the  wall;however,  our  approach  predicts  that 
the  layer  of  spurt  points  spreads  into  the  interior  of  the  channel  on  continued  loading. 
Layer  thickness  is  predicted  to  grow  continuously  in  loading  to  a  thickness  that  should  be 
observable,  provided  secondary  (two-dimensional)  instabilities  do  not  develop. 

Our  analysis  suggests  other  ways  in  which  experiments  might  be  devised  to  verify 
the  dependence  of  spurt  on  material  properties:  (i)  produce  multiple  kinks  with  spurt 
layer  separated  from  the  wall,  (ii)  produce  hysteresis  in  flow  reversal  (Fig.  9).  Our  model 
predicts  circumstances  under  which  a  different  path  can  be  followed  in  sudden  reversal  of 
the  flow  than  would  be  followed  by  a  sequence  of  solutions  in  which  the  pressure  gradient 
is  reduced  to  zero  and  reloaded  again  (with  the  opposite  sign)  to  a  value  of  somewhat 
smaller  magnitude.  Such  behavior  does  not  seem  likely  to  be  explainable  by  a  wail  effect. 

An  observation  that  has  been  associated  with  spurt  and  related  phenomena  is  the 
variation  in  kinematic  quantities  that  appears  to  be  periodic;  such  variations  are  inferred 
from  the  appearance  of  extrudates  at  the  exit  of  a  die  or  capillary  in  which  spurt  (or  slip,  or 
other  wall  effect)  is  presumed  to  be  occurring.  The  analysis  in  Sec.  6  does  not  predict  any 
truly  periodic  behavior,  but  the  spiraling  into  the  spurt  attractor  is  shown  to  have  a  very 
high  frequency  and  slow  decay  rate,  so  that  for  practical  polymer  systems,  it  may  very  well 
appear  to  be  periodic  and  persistent.  Furthermore,  since  these  oscillations  involve  large 
changes  in  the  normal  stresses,  they  could  very  well  lead  to  “periodic  distortions”  in  the 
stream  extruded  from  a  die  in  which  such  behavior  is  occurring.  Our  prediction  is  that 
the  average  “period”  of  such  distortions  would  be  proportional  to  the  viscosity  ratio,  e; 
however,  because  of  the  —  <r / e  term  in  the  denominator  of  Eq.  (6.14),  the  apparent  period 
could  vary  significantly  in  time  until  the  magnitude  of  the  a  contribution  to  the  shear  stress 
is  small  compared  to  s.  The  “ultimate  period”  is  predicted  to  be  inversely  proportional  to 
T,  to  leading  order. 
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Even  if  it  can  be  established  that  spurt  is  a  distinct  phenomenon  from  wall  slip  (or 
other  wall  effects),  our  model  is  not  the  only  possible  alternative.  Others  have  proposed 
that  spurt  is  related  to  non-monotone  constitutive  behavior  (see  for  example  Refs.  [6],  [15], 
or  [12]).  The  features  of  hysteresis  deduced  in  Sec.  5D  stand  in  marked  contrast  to  other 
plausible  predictions  of  the  nature  of  the  hysteresis  in  spurt  [12];  experiments  could  be 
devised  to  verify  which  theory  is  correct. 

Finally,  the  most  important  and  perhaps  the  easiest  experiment  to  perform  to  verify 
our  theory'  is  to  produce  latency.  Eq.  (6.4)  predicts  long  latency  times  for  data  correspond¬ 
ing  to  realistic  material  data;  no  sophisticated  timing  device  would  be  required,  nor  would 
the  onset  of  the  instability  be  hard  to  identify.  The  increase  in  throughput  is  predicted  to 
be  so  dramatic  that  simple  visual  inspection  of  the  exit  flow  would  probably  be  sufficient. 

8.  Conclusions 

Phase-plane  and  small  parameter  asymptotic  techniques  have  been  used  to  analyze 
single  and  two  relaxation-time  models  that  approximate  the  behavior  of  highly  elastic 
polymeric  liquids  in  shear  flow  through  slit  dies  driven  by  a  pressure  gradient.  Results  of 
the  analysis  show  that  both  models  exhibit  the  same  qualitative  features  for  realistic  data 
based  on  polymer  melts.  The  results  provide  an  explanation  of  the  spurt  phenomenon 
observed  experimentally;  this  explanation  differs  from  “wall  slip”  in  that  it  stems  solely 
from  material  properties  of  the  polymer  and  not  from  the  interaction  of  the  polymer  with 
the  wails.  In  addition,  the  analysis  predicts  new  phenomena  that  should  be  observable 
in  new  experiments:  latency,  shape  memory,  and  hysteresis  in  cyclic  loading  and  loading 
that  induces  sudden  flow  reversals.  The  key  to  understanding  the  dynamics  of  the  approx¬ 
imating  systems  is  fixing  the  location  of  the  discontinuity  in  the  strain  rate  induced  by  the 
non-monotone  character  of  the  steady  shear  stress  vs.  strain  rate. 

Although  our  analysis  applies  only  to  the  special  constitutive  models  we  have  studied, 
we  expect  that  the  qualitative  features  of  our  results  appear  in  a  broad  class  of  non- 
Newtonian  fluids.  Our  analysis  has  identified  certain  universal  mathematical  features  in 
the  shear  flow  of  viscoelastic  fluids  described  by  differential  constitutive  relations  that  give 
rise  to  spurt  and  related  phenomena.  Since  we  have  seen  that  the  qualitative  features  of 
two  models.  ( J SO ),  and  ( JSO2 )  are  very  similar  in  the  range  of  interest,  we  choose  to 
describe  the  universal  features  by  abstracting  ( JSO ).  The  key  feature  is  that  there  are 
three  widely  separated  time  scales,  each  associated  with  an  important  non-dimensional 
number  (a.  s,  and  1,  respectively),  when  scaled  by  the  dominant  relaxation  time,  A-1. 
Each  of  these  time  scales  can  be  associated  with  a  particular  equation  in  system  (JSO): 

rparab  :=  ph2X/ fJ.  &  —  rparab^ 

Tshear  •’=  Tj/fi  £  =  rshe^T\  (8.1) 

^normal  ;=  A  1  1  =  ^normal  • 

Note  that  r,hear  can  be  associated  with  with  £i2/^2>  a  secondary  shear  viscosity  in  sys¬ 
tem  ( JSOz  )•  We  assume  that  rparab  <  rshear  <  rnormai,  as  in  experiments  of  Vinogradov 
et  a 1.  [17]. 

To  see  how  these  time  scales  are  associated  with  specific  equations,  it  is  convenient  to 
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rewrite  system  ( JSO)in  the  following  equivalent  form: 


avt  -  <rx  =  £vxz  +  / 

T  +  a/j1  vt  ( x',t)dx ' 
£ 


ct  —  (Z  -h  1) 


(Z  +  1) 


+  1 


T-a  +  afzvt(x',t)dx' 

+  cr— - - -  =  —  Z  . 

£ 


(8.2) 


Setting  a  =  0  in  system  (8.2)  yields  the  approximating  system  (3.1):  the  extreme  smallness 
of  T-parab  allows  us  to  disregard  inertial  effects  of  the  first  equation.  So  long  as  Z  +  1  is  not 
0(c),  therefore,  the  effective  relaxation  time  in  the  second  equation  is  e.  Consequently,  cr 
achieves  its  steady  or  latent  value  on  a  time  scale  of  order  rs hear.  Finally,  as  a  consequence 
of  the  third  equation,  Z  tends  either  to  a  classical  steady  state  value  (less  than  —1/2),  or, 
in  case  of  spurt,  to  a  value  near  —1  with  effective  relaxation  time  1.  We  conclude  that  it  is 
this  general  structure  of  evolution  equations  governing  the  flow  that  causes  the  approach 
to  spurt  to  have  three  distinct  phases. 
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